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We numerically study a billiard system with a time-dependent force, and our results suggest 
the existence of a limitation on possible transitions between steady states in Hamiltonian chaos, 
in analogy to the limitation on transitions between equilibrium states described by the second law 
of thermodynamics. This limitation is expressed in terms of irreversible information loss, which 
is defined for each trajectory through Lyapunov analysis. As a key step in the study, we demon- 
strate a non-equilibrium equality which means that the average of the inverse exponential of the 
irreversible information loss is unity, where the average is taken over initial conditions sampled from 
. the microcanonical ensemble. 

(N 

£ ■ 05.45.+b,05.70.Ln 



3 



The second law of thermodynamics is a formalization of a fundamental limitation on possible transitions among 
equilibrium states . Although the validity of the second law has been confirmed conclusively without exception for 
£>^j \ more than a century, there yet exists no clear understanding of how this limitation emerges from purely Hamiltonian 
systems. 

When the dynamical system in question possesses a mixing property, time correlations of dynamical variables decay. 
This decay characterizes the directional evolution in dynamical systems. For this reason, it may be natural to believe 
U: that the chaotic behavior of a system is related to the limitation on possible transitions between its steady states. 
<—| ■ However, as far as we know, such a description has never been formulated explicitly. 

In this Letter, we study non-steady behavior of Hamiltonian chaos resulting from a change in the value of a 
parameter during a finite time interval. We define a quantity /, called the 'irreversible information loss', for each 
trajectory, and demonstrate that it satisfies the non-equilibrium equality 
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(exp(-7)) mc - 1, (1) 



where (-} mc denotes the average over initial conditions sampled from the microcanonical ensemble on an energy surface 
. specified initially. Through the Jensen inequality, (exp(— I)) mc > cxp(— (I) nlc ), the inequality 

W mc >0 (2) 



is derived. This may be regarded as a representation of the limitation on possible transitions. 

Our search for a non-equilibrium equality of the form ([!]) for Hamiltonian chaos was initially motivated by its 
obvious relation to the Jarzynski equality 

« : (exp(-/3(W - AF))) c = 1, (3) 

> 



where W is the work done by an external agent, AF is the free energy increment for a state transition, and (•) is the 
average over all possible histories, each of whose weight is determined by the probabilities of equilibrium fluctuations 
of the system in an isothermal environment with inverse temperature (3. We note that the minimum work principle 
(W) c — can be derived as the result of (||). Considering the similarity of (Q) and it is natural to call ([!]) 'a 
Jarzynski- type equality in Hamiltonian chaos'. 

Although the analysis we give here can be applied to a fairly general class of Hamiltonian systems, we focus on a 
billiard system with a time-dependent external force. Here, each phase space point T is specified by the canonical 
coordinates (<7i,<Z2jP1jP2)j and the Hamiltonian we study is given by 

H({quPih /) = ~(p? +Pl) + \{r(<li) 2 + ql) d/2 + /<?!, (4) 

with 

r(qi) = 6>(<2i)max(<7i - a, 0) + 0(-gi)min(gi + a, 0), (5) 

where 9 is the Heaviside step function. The quantity / in (^) is an external force and the potential represents a soft 
interaction with a stadium-shape wall. In this Letter, we set (d,a) = (8.0,0.5). We numerically solved the evolution 
equations 
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dqi , , 

~dt = (6) 
dpi = _dH_ 

dt dqi 1 ' 

by employing the 4-th order symplectic integrator method with a time step St = 1CP 3 . 

We begin our discussion by describing the basic properties of the system in the case of a time-independent force, 
say / = /o- Then, the energy is a constant, say Eq. As an example, we consider the case in which Eq = 1.0 and 
fo = 0.0. Our numerical results for this system lead us to conclude that it does indeed possesses the mixing property 
with respect to the microcanonical measure. We thus assume that a continuous distribution of initial conditions on 
an energy surface evolves into the microcanonical ensemble in the weak sense when the system is left for a sufficiently 
long time. 

The chaotic nature of the system is explored by applying a linear analysis (the 'Lyapunov analysis') to trajectories. 
Through the numerical integration of the linearized equations of motion corresponding to (^J) and (Q), 

~dT = Spi ' (8) 

we obtained a tangent evolution map from time s to t, T(t, s). 

To make the discussion more concrete, we choose a specific reference time and consider a set of four orthogonal 
unit vectors {ei(t*)} defined in the tangent space at the phase space point r(t*) at time t*. Then, the vector ei(t*) 
evolves to T(t, t*)ej(£*) in the tangent space at the phase space point T(t). By employing the Gram-Schmidt procedure 
T(i,i*)ej(t*) (where t > t„) can be expressed as 

4 

T{t,U)e t {U) =5^e 3 -(t)L Ji (t,t.), (10) 
i=i 

where Lij is the element of an upper triangular matrix whose diagonal elements are positive, and {ei(t)} is the 
set of four orthogonal unit vectors defined in the tangent space at the phase space point T(t). The ith Lyapunov 
exponent A, is calculated as the long time average of the quantity A^ (t) given by 

\i(t) = - log L u (t,U). (11) 
at 

Although the vector e;(i) and the quantity Xi(t) depend on the choice of the set {ei(t*)}, the value of Ai is independent 
of this choice for the system in question, and we have Ai = — A4 > and A2 = A3 = 0. 
There is a special vector e^ S ' ) (r(t); fa) at the point T(t), which is defined as 

e< S) (r(t);/ )= lim ei(t). (12) 

t*— > — OO 

In our numerical experiments, we insure that e\ (T(t): fo) is represented by ei(t) by choosing the value of t* so that 
t — t* is sufficiently large for the condition 

|l-( ei (*),ei(*)) 2 |<e (13) 

to be satisfied for two vectors ei(f*) and e'^t*) chosen randomly, where e is arbitrarily chosen to be 10~ 6 . The vector 
e^\r(t); fo) corresponds to the most unstable direction and is called the 'first Lyapunov vector' at the point T(t). 
Using this vector, we can define the first local expansion ratio as 

xf ) (T(t);f ) = ^\ogLf 1 \t,U), (14) 

(S) (S) 

where L lx (t,s) is the magnitude of the vector T(t,t*)e\ (r(t*); fo). 

One may expect that e^\T(t); fo) and A| (r(i);/o) can be defined similarly to e^ S ' ) (L(t); fo) and (T(t); fo). 
However, numerically, the convergence condition in this case, taking a form similar to (131), is difficult to realize. 
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We thus added an artificial, small damping term —r/Spi to (0) in order to obtain e^\T{t); fo). We expect the value 

(S) 

obtained for e 2 (T(t); fo) in this manner to be well-defined in the limit 77 — > 0+ after first taking t* — > —00. We call 

(S) (S) 

the subspace spanned by e\ (T(t); f ) and e 2 (T(t); fo) the 'semi-unstable space'. 

Next, we consider the system with a time-dependent force. We consider a function f(t) that changes from fo to f\ 
only during a finite time interval [to,ii]. The initial conditions at t = r_ (< to) are sampled from the microcanonical 
ensemble on an energy surface with energy Eq. 

The Lyapunov analysis for the time-dependent system is developed in the following way. First, using the Gram- 
Schmidt decomposition of the linearized evolution map T(t, s), the time evolution of a set of orthogonal unit vectors 
is determined. However, the Lyapunov vector e\ (T(t); fo) does not have meaning for t > to- Instead, we define the 
'forward Lyapunov vector' e[ F \t) as the unit vector at time t obtained from ep' ) (r(r_); fo) under the Gram-Schmidt 

CB^ (SI 

procedure, and the 'backward Lyapunov vector' e\ (t) as the unit vector from which e\ (T(t + ); fi) is obtained under 
the Gram-Schmidt procedure, where t\ < r + . Since Li\J\t, r_) and L[f\r+,t) are determined simultaneously with 
ef^(t) and ej B ^(t) through these Gram-Schmidt procedures, we define the forward local expansion ratio X {F \t) as 

Af)(t) = |log4 F )(i,r_), (15) 

and the backward local expansion ratio x\ B ^ (t) as 

\T\t) = j t \ogL ( S\r + ,t). (16) 

When the force / is time-independent (say fo), three quantities x[ F \t), x[ B \t) and X^ (T(t); fo) are identical. This 
leads us to conjecture that the difference between A,- (t) and A ■ (t) provides a useful characterization of phenomena 
exhibited uniquely by time-dependent systems. Also, the long time average of A^ (F(t);/ ) is known to be equal to 
the information loss rate. With these in mind, we define the irreversible information loss / as 

1= lim lim [ T+ dt^lX^^-XfHt)}, (17) 

~ i— 1 

where we remark that the sum is taken over indices corresponding to the semi-unstable space, not unstable space. 

Now, we present numerical results which indicate that the non-equilibrium equality (|l]) holds for the irreversible 
information loss / defined by (|l7|). We first note that e| B (i) cannot be calculated directly, because the backward 
evolution of the semi-unstable space is not numerically stable. Thus, in order to be able to apply the type of analysis 
we discuss here to the general situation, we need to devise some more sophisticated algorithm to calculate e\ (t). In 
this Letter, however, we focus on the simplest case in which / is changed from fo to f\ instantaneously at t = to- 

(B) (H) (S) (S) 

In this case, when t > t , e\ (t) and A- Jt) are equal to e- ^(r(i);/i) and A- (r(i);/i), respectively. That is, the 
integration over the region to < t < t + in ( |17|) can be done numerically. The other part of integration can be evaluated 
with the aid of the time-reversed trajectory {T(t)} t ^Z T+ , which is the solution to the Hamiltonian equation with the 
time-reversed external force, {/(f)}t=_ r+ , and with the initial condition F(— t + ) = 1ZT(t + ), where 1Z is the operator 
that changes the sign of the momentum. We can carry out Lyapunov analysis for the time-reversed trajectory, and we 
express all quantities obtained in this analysis by simply adding a tilde to the corresponding quantities in the normal 
Lyapunov analysis. From the reversibility of the Hamiltonian equations, the integration over the region r_ < t < to 



in (17) turns out to be equal to 



where sin^^io) is the four dimensional volume spanned by the vectors e*f\to) and TZe\ F \to) {i — 1,2), and 
sin^( B ^(to) is defined similarly. The quantity ( [i"8| ) can be calculated by using the forward evolution of the semi- 
unstable space, because ^'(t) and (r B \t) are identical to A- S '(f (— t); f Q ) and e- S '(f (— t); f ). 

In this way, we can calculate the numerical value of / for each trajectory and check whether or not (exp(— J)) mc 
is unity. However, there is a slight complication here; trajectories with large negative / are rare, but they contribute 
greatly to (exp(— I)) mc . Thus, in order to avoid large numerical error in the computation of (exp(— 7)) mc that may 
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result from the inclusion of too many such rare trajectories, we divided the range of values of I into discrete intervals, 
and removed from consideration all trajectories whose values of I fell inside such intervals that contained fewer than 
a certain minimum number of data points. In this way, we obtained the normalized distribution function n(7), which 
is plotted in Fig. [I]. Using the data contained herein, for example, we obtain (exp(— 7)} mc = 1.03. 

However, the value obtained for (exp(— 7)) mc in this manner depends on how we divide the range of values of I 
into intervals. (In fact, without changing the shape of 11(7) shown in Fig. [I] significantly, we can make (exp(— I)) mc 
identically equal to 1.) It is thus necessary to find a better method to determine the value of (exp(— I)) mc from our 
data. One such method is to consider the asymmetry around the peak value of the function 11(7). As seen in Fig. [^, 
it is expected that the equality 

io g n(7)-io g n(-/) = / (19) 

holds. That is, II(I) should possess the symmetry described by the fluctuation theorem 0-0]. Then, since (0) can be 
derived from (|l9|), and since the situation depicted in fig. ^| seems to be independent of the manner in which we treat 
the data, we conclude that (Q) holds. 

Finally, we present a theoretical argument for ([!]). We note, however, that this argument lacks mathematical rigor. 

Let us consider M trajectories whose initial conditions at t = r_ are sampled from the microcanonical ensemble, 
where M is a sufficiently large number. Here, in order to simplify the theoretical arguments, we consider the micro- 
canonical ensemble to be defined by the uniform measure on the region between two energy surfaces with energies Eq 
and E + SE, where SE < E . 

We divide the phase space into small cells {A^}^, each of which has volume e. We also discretize time as 
t n = n(r + — T-)/N + r_, where < n < N. We can characterize a trajectory by a set of integers («o, • • • , in), where 
the phase space point on the trajectory at time t n is included in the cell Ai n . We then define the path probability 
P(io, • • • , «jv) as the ratio of trajectories characterized by (io, • • • , in) on the set of M. 

The probability p(i n , t n ) of finding the phase space point of an arbitrary one of the M trajectories in the cell A^ n 
at time t n , is given the sum of P(io, *iv) over all configurations (ii, ■ • ■ , i n _i, i n +i, iiv)- We note that p(i n , t n ) 
can be interpreted as the time evolution of the distribution function for M phase space points which are sampled 
from the microcanonical ensemble at t = r_ . 

When r + — t_ is sufficiently large, phase space points on the M trajectories at time may spread over the semi- 
unstable space in the ijyth cell, subject to the relation p(iN, tjsr) ^= 1. Therefore, it is expected that P(io, • • • , in) can 
be expressed as 

P(t , = A (B) (io,---,iN)p(iN,t N ), (20) 

where A^ B ^(zo, • • • ,ijv) is the probability of finding an arbitrary one of the trajectories that are in the semi- unstable 
space in the «Arth cell at time t = t/v to be characterized by («o, • • • ,in)- Then we introduce the similarly defined 
quantity A( F )(io, ■ ■ ■ ,%n) as the probability of finding an arbitrary one of the trajectories that are in the semi- unstable 
space in the ioth cell at time t — to to be characterized by (iq, ■ ■ ■ , in)- Then, using (pp|), we have the trivial identity 

,. ^. / ""'AW(i ,---,iN)p(iN,tN) 

Now, we conjecture that there is an appropriate limit in which A^ F ^(io, • ■ • , «tv) and A^ B ) (iq, ■ ■ ■ , z^r) approach 
exp(— Yl1=i (*)) an< l ex P(~ J T+ S?=i \' B ' ) (0); respectively. Also, from Liouville's theorem, we know that 
p{ioi to) = p(iNi iiv) holds in the limit e — > 0. Therefore, it is reasonable to expect that (|l|) can be derived from (pi]). 

In summary, we have presented numerical evidence for the validity of the non-equilibrium equality (|l|) together 
with a theoretical argument. 

In a previous paper, by numerically studying a Fermi-Pasta-Ulam model with a time-dependent nonlinear term 
we have found that the Boltzmann entropy difference has a certain relation to the excess information loss H cx . In the 
billiard model we study presently, H cx is given by the time integration of (t) ~ \\ (r(t); /i) over the region t > t - 
Therefore, I is expected to be related to the Boltzmann entropy difference in the thermodynamic limit. Indeed, it has 
been conjectured that these quantities are equal in the thermodynamic limit ||. Direct numerical evidence supporting 
this conjecture will be presented elsewhere, together with a complete theory. 

Our study is apparently related to the fluctuation theorem, as indicated by the fact that our data apparently satisfy 
(p"9|). However, it should be pointed out that the entropy production ratio in thermostatted models is given by the 
phase space volume contraction ratio, which is related to A^ F ^ — \\ F \ not A| F ^ — A^ B ^. Therefore, in extending our 



p,, ■ x A( F )(i ,--- ,i N )p{i ,t ) 

P (*o j ---,^)t7rT7^ —TT- 7 \ = 1- ( 21 ) 
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analysis to non-equilibrium steady states, the irreversible information loss defined by (|17|) does not become the entropy 
production. Rather, we believe that this comes to represent a certain quantity related to steady state thermodynamics, 
which was proposed phenomenologically by Oono and Paniconi |l0[ . 
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FIG. 1. The normalized distribution function 11(7). The values of Eo, fo, fi, r) and r± — to are here 1.0, 0, 0.5, 0.02 and 
±400, respectively. We construct a histogram by dividing the interval [—5, 5] into 50 boxes for the 9964 samples considered 
and removing boxes in which the number of data points is less than 25. 
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